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Method and arrangement for medical X-ray imaging_ 
BACKGROUND OF THE INVENTION 

5 TTiree-dimensional X-ray Imaging Is based on taking several 1-D or 2-D 
projection Images of a 3-D body from different directions. If l-D projection 
images are available fix)m all around a 2-D slice of the body with dense 
angular sampling, the inner structure of the slice can be determined. This is 
known as Cbmputerized Tomography (CT) imaging technology; which is 
10 widely used in medicine today. A crudal part of CT technology is the 

reconstruction algorithm taking the X-ray images as aigument and producing 
a voxel representation of the 3-D body. 

In many practical cases X-ray projection images are available only from a 
15 limited angle of view. A collection of X-ray images of a 3-D body is callecf 
sparse projection data if (a) the images are taken from a limited angle of 
view or (b) there are only a small number of images. Sparse projection data 
does not contain sufficient information to completely describe the 3-D body. 

20 However, some a priori information about the body is typically available 
without X-ray imaging. Combining this Information with sparse projection 
data enables more reliable 3-D reconstruction than is possible by using only 
the projection data. 

25 Traditional reconstruction algorithms such as filtered backprojecHon (I=BP), 
Fourier reconstruction (FR) or algebraic reconsfrucHon technique (ART) do 
not give satisfactory reconstructions from sparse projection data. Reasons for 
this include requirement for dense full-angle sampling of data, difficulty to 
use a priori information for example nonnegativlty of the X-ray attenuation 

30 coefficient and poor robustness against measurement noise. For example the 
FBP method relies on summing up noise elements with fine sampling, leading 
to unnecessarily high radiation dose. 
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BRIEF DESCRIPTION OF THE INVENTION 



The aim of the invention is to overrome the problems met in 3-D 
5 reconstruction of the body that occur when using traditional reconstruction 
algorithms with sparse projection data. This Is achieved by a method for 
producing three-dimensional Infbmiation of an object In medical X-ray 
Imaging Ip which method the object is modelled mathematically 
Independently of X-ray imaging. The object Is X-radiated from at least two 
10 difterent directions and the said X-radiation is detected to form projection 
data of the object Said projection data and said mathematical modelling of 
the object are utilized in Bayesian Inversion based on Bayes' formula 



p(x\m) = 



p(m) 

to produce three-dimensional information of the object, the prior distribution 
15 pp,<x) representing mathematical modelling of the object, x representing the 
object Image vector, which comprise values of the X-ray attenuation 
coeffident Inside the object, m representing projection data, the lilcelihood 
distribution p(m|x) representing the X-radiation attenuation model between 
the object image vector x and projection data m, p(m) being a normalization 
20 constant and ttie posteriori distribution p(xim) r^resenting the three- 
dimensional information of the object. 

The invention also relates to a medical X-ray device arrangement for 
producing tiiree-dimensionai information of an object in a medical X-ray 
25 Imaging, said medical X-ray device arrangement comprises: 

- means Ibr modelling the object mathematically independentiy of X-ray 
imaging 

- an X-ray source for X-radiating the object from at least bNo different 
directions, 

30 - a detector for detecting the X-radiation to form projection data of tfie 
object. 
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- and means for utilizing said projection data and said mathematical 
modelling of the object in Bayesian inversion based on Bayes* formula 



Ppr i^)p{rn I X) 



to produce three-dimensional infomiation of the object, the prior distribution 
5 ppr(x) representing mathematical modelling of the object, x representing the 
. . object Image vector, which comprises values of the X-ray attenuation 
coefficient inside the object, m representing projection data, the iikeKhood 
distribution p(m|x) representing the X-radlation attenuation model between 
the object Image vector x and projection data m, p(m) being a normalization 
10 constant and the posteriori disbibution p(x|m) representing the three- 
dimensional information of the object 

The invention Is based on that biological issues have that kind of statistical a 
priori information that this information can be utili^ succesfuliy vwth 

15 bayesian Inversion in medical x-ray Imaging, The suitable a priori Infomiation 
makes possible to model the biological tissue mathematically accurately 
enough and Independentiy of X-ray Imaging. From biological tissue it is 
possible to compile kvalltative stiijcturai information which makes It possible 
to utilize tiie bayesian meOiod successfully to solve the problems In medical 

20 Hiree-dimensional x-ray imaging. There are certain regularity In biological 
tissues and this regularity Is useful especially witii ttie bayesian method. 

Fbr example 10 x-ray images Is taken from breasts of different persons. 
From these x-ray Images is noticed that there Is much similarity In the 
25 statistical stiucture of the breasts between different people. In other wonjs 
biological tissues and x-ray images taken from Oie biological tissues has 
similar or almost similar statistical structure between different persons. 



30 



Witii Bayesian inversion it is possible to utilize a priori information efficientiy 
in 3-D reconstruction from sparse projection data, because the suitable a 
priori information fi-om tiie biological tissue makes possible to model the 
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biological tissue mathematically aocurately enough and independently of X- 
ray imaging. Any collection of projectfon data can be used in the 
reconstruction. Application-dependent a priori knowledge can be used to 
regularize the ili-posied reconstnjction problem. 

5 

TTiis Invention Improves the quality of 3-D reconstrucHqns over traditional 
methods. In addition, the number of radiographs can be minimized without 
compromising quality of the recx)nsbiictlon. This is very important in medical 
applications since the X-ray dose of the patient can be lowered. 

10 

BRIEF DESCRIPTION OF THE HGURES 

Figure 1 shows a simple example of X-ray Imaging. 

15 Figures 2A -2B show a difference between global tomography and local 
tomography. 

Figures 3A - 38 show different types of sparse projection data. Every black 
dot represents a locabon of the X-ray source for taking one 
20 projection image. 

Rgures 4A - 4B show examples of parts of boundaries and cracks of an 
object visible and invisible without a priori informatton. 

25 RgureS illustrates "pendl beam" X-ray attenuation model. 

Rgure 6 Is basic flow chart of the method according to the Invention. 

Rgure? Is an intraoral X-ray device arrangement presenting one 
30 preferred embodiment of the invention. 
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Rgure 8 shows a measurement geometry of dental llmrted-angle 
tomography with a digital intraoral sensor. 

Rgure 9 shows results of 3D reconstruction of head phantom In the first 
5 preferred embodiment 

DET/ULED DESOUPTION OF THE INVENTION 

In practical Imaging situations X-ray imag^ are not always available from all 
10 around the body. The body might be visible only from certain directions due 
to imaging geometry. For example tills is the case in 3-D mammography with 
ttie breast compressed against the detector, or In Intraoral dental imaging 
with the detector inside the patient's mouth. This situation Is called limited- 
angle tomography. Also, the region of Interest might be surrounded by tissue 
15 that need not be Imaged, like In extraoral dental imaging. This situation is 
called local tomography. In addition tiie number of radiographs should be 
minimized In medical applications for reducing the X-ray dose of the patient. 

In the preferred embodiments of tiie invention Bayesian inversion algorithms 
20 are used to create a new type of 3-D medical X-ray Imaging. It is 
Intermediate between a projection radiograph and a full computed 
tomography scan. Two steps are needed to perfomi successfully such 
Imaging: In step one, tiie doctor (a) collects all tiie prior information he or 
she has on tiie tissue of interest, and (b) takes tiie minimum number of 
25 radiographs containing Oie desired 3D information. In step two, a 
tomographic algorltfim based on Bayesian inversion Is used for 3D 
reconstruction. The algorithm contains bwo separate matiiematical models. 
Rrst, all a priori information (I.e. Infomiation tiiat is independent of X-ray 
imaging) is used to model the unknown object maOiematically. The model is 
30 put into tiie form of prior disbibution in Bayes' fonnula. Second, the 

measurement is modelled matiiemafa-cally. This Involves tfie geometry of tiie 
Imaging device, position of tiie detector and X-ray source during the 



wo 2004/019782 PCT/FI2003/0(H)625 

6 

exposure of each projection image and a model for X-ray attenuation inside 
dssue. TTie mathematical model of the measurement is put into the fomn of 
likelihood distribution in Bayes' formula. 



5 In figure 1 is shown a simple example of X-ray imaging, where an X-ray 
source 2 is placed on one side of an object 4 under Imaging. Radiation 
passes through the object and Is detected by a detector 6 on the other side. 
The X-ray source is for example a X-ray source of an intraoral X-ray source 
of a dentist of a dental panoramic X-ray device, of a surgical C-ami X-ray 
10 device, of a mammography device or of any other medical X-ray device and 
the detecbsr 6 is a detector of some of those devices. Usually the detector 6 
Is a digital isensor that can be thought of as a 2-D array of almost pointlil<e 
detectors. 

15 The 3-D body under imaging is modelled by nonnegative X-ray attenuation 
coefficient The value gives the relative intensity loss of the X-ray travelling 
wittiin a small distance dr: 

^ = —x(r)dr 

(1) 

20 The X-radiation has initial intensity when entering the object 4 and a 

smaller intensity when exiting the object The following equation shows 
the attenuation law: 

f f I'M 

/ x(r)dr = - / -jTrdr = log Jo - log Ji 



(2) 
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Where initial intensity Is known by calibration, and Intensity after object 

^1 Is known from the corresponding point value in a projection image. Thus 
the measured data is the Integral of x along the line L 

5 In ttie above model it Is not taken Into account (a) scattering phenomena 
resulting In X-rays changing their direction, or (b) the dependency of 
attenuation on the X-ray spectrum resulting In low-energy photons being 
more easily attenuated than high-energy ones. Effect (b) causes errors in 
measurements and Is sometimes referred to as beam hardening. More 
10 detailed models can be utilized in Bayeslan inversion to. overcome (a) and 
(b). 

In medical imaging the geometrical ai-rangements of the X-ray source and 
digital sensor vary according to the diagnostic task and equipment Rgures 
15 2A-2B and 3 illustrate the types of tomographic data resulting from different 
imaging situations. For darity, here are presented two-dimensional 
examples; similar situations can be considered in 3-D. 

In figure 2A-2B are shown two cases according to whether the whole object 
20 4 is fully visible In each projection or not These cases are called global 

tomography and local tomography, respectively. In figure 28 is radiated only 
ROI (Region Of Interest) 8. 

The types of data described In figures 3A - 38 cover a large range of specific 
25 imaging tasks. The choice of data collection dictates what kind of features 

and dietails It Is easiest to reconstruct reliably from the data. Thus the choice 

of projection views must be made according to the diagnostic task at hand. 

in traditional CT Imaging, projections are taken from all around the object 4. 

In the preferred embodiments of the Invention radiation dose is lowered by 
30 sampling the angular variable more sparsely. Figure 3A presents more dense 

angular data witii full angle and limited angle measurements. Rgure 38 
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presents more sparse angular data, also with full and limited angles of 
' measurement 



In figures 4A and 4B there are presented etamples of parts of boundary and 
5 cracks that are visible or indetedable in reconstruction without a priori 
information. Rgure 4A there is an object 4 under imaging with edge on the 
surface of the object The edge Is detectable on itie leftmost position where 
the edge is more parallel to the direction of the X-rays. In the figure 4B there 
is a cradc itiside of the object 4. The crack is detectable in the leftmost 
10 position parallel to the angle of mea^rement 

The main idea in Bayesian inversion method is to consider the Inverse 
problem as a problem of statistical inference. All variables are redefined to be 
random variables. The randomness reflects uncertainty of their actual values 
15 and the degree of uncertainty is coded In the probability distributions of 
these random variables. 

When modelling the measurement mathematically the X-ray attenuation 
model and the observation can be assumed linear. The model is presented in 
20 the following equation. 

(3) m=Ax+e, 

where the variables m, x and e are random variables. 
The matrix A comes from the pencil beam model for the X-ray attenuation. 
25 This is presented in figure 5. There the unknown 3-D body of the object 4 is 
divided into small voxels 10, and the matrix A contains tiie lengths of the 
path of the X-ray inside each voxel. In this way the integral In formula (2) is 
approximated with a simple numerical quadrature. 

30 Assume now that the object image vector x and the noise are jointiy 

independent random variables. TTie conditional probability distribution of x 
given the measurement m is given by Bayes' formula in the form 
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I ^ Ppria^)pnoise{m-Ax) 

p(m) 

(4) 

where p(m) is a normalization constant TTie density p(x|m) is called the 
posteriori distribution of x. Pnoise(m-Ax) = p (m|x) is a likelihood distribution 
representing the X-radiation attenuation model between the object image 
5 vector X and projection data m. ;>i; 

The density ^f" , called the prior distribution of x, is designed to contain all 
possible infomriation available of the actual object 4 Independently of X-ray 
Imaging. It is crudal, in contrast to several classical regularizatlbn methods, 
10 that the choice of the prior disbibution should not be based on the projection 
data. The rule of thumb in the design of prior distributions is that typical 
image vectors (say, of some existing library) should have high prior 
probability while atypical or impossible ones should have low or negligible 
probability. 

15 

In the framework of the Bayesian Inver^on theory, the posterior disb-ibution 
In formula (4) represents the complete solution of the 3D reconsbuction 
problem. 

20 To produce an Image of the ol)|ect 4 based on the posterior distribution, 
several alternatives ©dsL TTie most common ones are the maximum a 
posteriori estimator (MAP) and conditional mean estimator (CM). They are 
defined by the formulas 

25 P(^MAP I rn) — majcp{x | m), 



(5) 
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where the maximum on the right hand side is talcen over all x, and 
3^CM == j xp(x I m)dx. 

5 (6) 

Finding the MAP estimator is an optimization problem whlle finding the 01 
estimator is a problem of integration, 

In figure 6 is presented a basic flow chart of the method according to the 
10 invention. In method step 600 the object is modelled mathematically 

independently of X-ray imaging. In m^od step 602 ttie object is X-radlated 
from at least two different directions. In method step 604 the said X- 
radiation is detected to form projection data of the object In triethod step 
606 said projection data and said mathematical modelling of the objed: are 
15 utilized in Bayeslan inversion based on Bayes formula to produce three- 
dimensional information of tiie object. The meUiod step 600 Is also possible 
to perform after method step 602 or 604. 

In the first preferred embodiment of the invention is presented an application 
20 to dentel radiology. 

X-ray projection images are conventionally used In dental radiology. 
However, certain diagnostic taslcs require more precise knowledge of the 3D 
structure of tissue than is available in two-dimensional radiographs. Such 
25 tesks include implant planning and detection of bone loss between tooth 
roots. 

In figure 7 is presented an intraoral X-ray device 5 arrangement presenting 
the first preferred embodiment of the invention. It is important to notice that 
30 this is only an example of the medical X-ray device 5 arrangement where the 
Invention is possible to be utilized. 
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The medical x-ray device 5 In the preferred embodiments of the Invention Is 
for example a dental panoramic X-ray device, a surgical C-arm X-ray device 
or a mammography device. 



5 In the first preferred embodiment of the Invention the articulated arm 
arrangement 3 moves the X-ray source 2 to the right position. The X- 
radiation begins by pressing the exposure button 12. The X-ray source 2 X- 
radlates the object 4, which is for example teeth of a patient The detector 6 
detects ttie X-radiation. The Image Information whidi.is got by detecting the 

10 X-radiation is sent by communication link 16 to the computer 14. The 

computer comprises the software means 15 to process the image information 
acojrdlng to the invention. There can be more than one computer 14 and 
also the software means 15 can sitoate in more than one computer 14. For 
example the first computer 14 is computer which Is used In x-ray imaging. 

15 The second computer 14 is computer which Is used in processing the Image 
infonmation according to the invention. It is pos^ble to have the second 
computer 14 far away from the actoal medical x-ray device 5. I=tor simplicity 
in figure 7 is sliown only one computer 14. 

20 In the first preferred embodiment (Rg. 7) of the invention, the dentist's X- 
ray equipment Is used for teWng a set of 2D projection images that are used 
as input for Bayeslan 3D reconstruction algorithm. Such equipment includes 
an intraoral X-ray unit and a digital intraoral sensor. 

25 Benefits of this approach over conventional CT scan are 

- low cost and convenient usage, 

- high resolution of the projection Images, 

- possibility to teke as few radiographs as Is needed to capture the relevant 
3D information, minimizing radiation dose to the patient, 

30 - the possibility to choose Imaging directions so that the X-rays do not pass 
through the whole head but only the interesting tissue, furtfier reducing 
dose. 
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In the preferred embodiments of the invention, Bayesian inversion is used 
for the 3D reconstruction. Input data for the algorithm is the set of 
projection images and the following a priori knowledge: 
5 (a) Dental tissue consists of few appro}dmately homogeneous regions witti 
well-defined boundaries, 
(b) Dental tissue can only attenuate X-radiation, not intensify it 

The example of the detector 6 used In the first prefemed embodiment of the 
10 Invention Is based on charge coupled device (CCD) technology and has 
dynamic range of 4096 gray levels. TTie size of the acHve imaging area is 
34 mm * 26 mm and the resolution Is 872 * 664 pbcels. After exposure, each 
pixel contains an Integer which is proportional to the number of X-ray quanta 
tiiat hit the pixel's area. 

15 

Alternative detectors include any other digital intraoral sensor, digitized X-ray 
film, or any intraoral sensing device converting detected X-ray photons to a 
digital image. 

20 In the measurement geometry the focal point of the X-ray source 2 moves 
horizontally on a circle In the plane 10, with center at the detector 6, see 
Rgure 8. In Rgure 8, the detector 6 corresponds to one row of CCD pixels in 
the CCD sensor. 

In figure 8 is presented an example of taking seven radiographs of some 
25 teeth of a patient The total opening angle of the projection views is 55 
degrees in this example. Rgure 8 represents the xy plane, and the z- 
coordinate a)ds is thought to be perpendicular to the plane of the paper. 

In the preferred embodiments, mathematical modelling of the object 4 can 
30 be done with a prior distribution of the form in formula (7): 
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Ppr(a?) = exp(-a (»)) 

N 

(7) 

where tiie sum is taken over a collection of 3D neighbourhoods N and the 

value depends only on the values of voxels belonging to the 
neighborhood and a is a positive reguiarization parameter used to tune 
the width of the prior distribution. The neighborhood N is t/pically a set of 
voxels whose center points are closer to each other than some predeifined 
maximum distance. 



It Is to be noted that formula (7) does not define a probability disb-ibution 
since it does not integrate to 1. Instead, the integral diverges. However, 
when used in the Bayes formula (4) with a suitable lil<elihood disb-ibution, 
the posterior distribution is integrabie. The same remaric concems the priors 
15 derived from (7), that is, formulae (9), (10) and (11). 



In the preferred embodiments of the invention the three-dimensional 
problem can be reduced to a stack of two-dimensional problems each 
corresponding to a plane determined by a constant value of z. Here figure 8 
20 represents exactiy the situation in the xy plane, i.e. the plane z=0. Each row 
in the detector corresponds to one such 2D problem. This approach leads to 
small approximation error because the X-ray source does not move in the 
correct plane for nonzero z, and this error is neglected. 
Next the modeling of a 2D tomographic problem is explained. Let 

25 

(8) ma) = Axa) + ea) 

denote the jth 2D tomographic problem. Here the vector mQ) contains the 
readings on jth row from each of the seven radiographs. The vector x(j) Is 
30 the jth slice of the 3-D representation x of the object 4 under imaging. The 
matrix A is the pendl beam model matrix for the 2-D tomography problem. 
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Now xQ) is a 2D array of pixels. The pixels are denoted by x(j)[k,q]/ where 
k=l,2,...,K is row Index and q=l,2,...,Q is column Index. 

TTie mathematical modelling of the object 4, I.e. incorporation of prior 
5 infonnation, is next explained for the 2D slice. Define 

K Q 

p,r{^^^) = eKp(-a^J2U{N{x^%q]))) 

(9) 

10 where the value U(N(x(j)[k,q])) depends only on the values of pixels 

belonging to the neighborhood N of the p'Dcei xQ)[k,q], the constant c is a 
normalization constant and a is a positive regularlzation parameter used to 
tune the width of the prior distribution. The choice of the neighborhoods 
N(x(j)[k,q]) is arbitrary, but a typical choice is the set of pixels whose 

15 cenlBrpolnts are closer to the centerpoint of xO)[k,q] than some predefined 
distance. 

In the preferred embodiment, the neighborhoods are chosen to consist of 
two adjacent pixels. Further, the function U is diosen to be a power of 
20 absolute value of the difference between adjacent pixels. These choices lead 
to formula (10): 



25 Pp.(^^'^)= exp(-a(X;f:|a;(^>lfc,g]-^^''>[fc + l,9]r + 

+ EE - q + l]f )) 

(10) 

30 where s Is a positive real number, and a is a regularlzation parameter used 
to tune the width of tiie prior distribution. 
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In the equation (11) Is presented non-nonnalized total variaton (TV) 
distribution, wlien s=l in equation (10). This is the case according to the 
preferred embodiment of the invention, because when s=l the prior model 
becomes such that it gives high probability for objects consl^'ng of a few 
5 areas of different attenuation with well-defined boundaries. 

piTV(a:«)) = exp ( - a( X; E l^^'^f*' " ^^'^t^b + l,q]\ + 



10 



ifc=l q=l 



(11) 

X-radlatioh an only attenuate, not istrengthen. Inside tissue. This leads to 
posltivit/ prior pos defined by 

15 (12) pos(x(j)) = Iff all pixels of xQ) are positive, 0 otherwise. 

The explanation of modeling of one 2D tomographic problem is now 
complete. 

20 The stack of 2D problems are connected to each other by demanding that 
consecutive 2D slices x(j) and xQ-l) should not be very different 
Mathematically this Is expressed as follows: 



25 



(13) pr3D(xa)) = exp(- y SZ I xa)ik,q] - xa-l)[k,q3 |), 

where the sums are taken over k=l,...,K and q=l,...,Q and y > 0 is another 
regularization parameter. So the prior distribution in fonnula (4) Is the 
product of (11), (12) and (13): 



30 



(14) ;^'^(xa)) = prTV(xa))pos(xa)) pr3D(xa)). 
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The measurements are taken into account In the form of the likelihood 
distribution 

At 

5 (15) 

where 2 is the covariance matrix of the Gaussian noise vector e. cIs a 
normalization constant. 

in reality, X-radlatlon measurement noise Is Polsson distributed. Formula 15 
10 uses a Gai^lan approximation tx) the Poisson distribution. 

All the parts are presented of the right hand side of (4). Next it Is looked for 
the MAP estimate of xQ), that Is, the Image xO) tiiat gives the largest value 
for the posteriori distribution (4). This Is equivalent to finding the Image xQ) 
15 that minimizes the expresaon 

At 

- log(prTV(a;(^>)) - log(pr3D(aj^>)) 

(16) 

20 with the additional requirement that ever/ pbcel of x(j) Is positive. 

Minimization of (16) is difficuit since F Is not a differentiable function due to 
the absolute values present in (11) and (13) and due to the sharp cutiang in 
(12). The absolute values In (11) and (13) are replaced by 



25 



\t\ « h^{t) = ilog(cosh05t)), 



(17) 

where p>0> enabling the use of efficient gradient-based minimization 
30 routines. Tlie posltivlty constraint Is taken care of by an Interior search 
minimization method. The results are shown in figure 9 in comparison with 
traditional tomosynthesis technique based on backprojection. 
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In this description explained mathemaHc modelling and other image 
Information processing are perfomied by software means 15 to accomplish 
three-dimensional information of the object 4. 

5 

in the second preferred embodiment the Invention Is utilized In dental 
panoramic x-ray Imaging. The object 4 Is typically teeth of the patient and 
the medical x-ray device 5 anrangement Is dental panoramic x-ray device 
arrangement , . 

10 

In the third preferred embodiment the Invention is utilized In mammography, 
for example in l=FDM (Full Reld Digital Mammography). TTiere the object 4 Is 
breast of human being and the medical x-ray device 5 arrangement is 
mammography de>^ce arrangement 

15 

In the fourth preferred embodiment the invention is utilized in surgical x-ray 
imaging. There the object 4 is the body part under surgery and the medical 
x-ray device 5 arrangement is surgical x-ray device arrangement 

20 Jn the first, second, third and fourth preferred embodiments of the Invention 
the basic metiiod steps are same as mentioned with the flow chart in figure 
6. The utilizing of the invention in the second, third and fourth preferred 
embodiment Is similar to what Is described with the first preferred 
embodiment of the invention and elsewhere in this application except 

25 different medical x-ray imaging applications and their differences because of 
different medical x-ray devices and different objects to be x-ray Imaged. 

Although the Invention is described above with reference to the examples 
illustrated In the attached figures, it Is to be understood that the Invention Is 
30 not limited thereto but can be varied In many ways within the Inventive idea 
disclosed in the attached dalms. 



